# REPLICATION CODE FOR "POLITICAL DISTINCTIVENESS AND DIVERSITY AMONG LGBT AMERICANS"
# PHIL JONES
# 9.30.20

library(car)
library(foreign)
library(ggplot2)
library(survey)

set.seed(2020)


######################################################################
# ANALYSIS OF JOURNAL ARTICLES (TABLE 1)
######################################################################

articles <- read.csv("journal_articles.csv")
round(colMeans(articles[,c(7:12)], na.rm=T)*100,1)
summary(articles)



######################################################################
# LOAD AND WEIGHT SURVEYS, GIVE NAMES OF DVS
######################################################################

load("anes.RData")
w.anes <- svydesign(~1,data=anes,weights=anes$weight)
load("cces.RData")
w.cces <- svydesign(~1,data=cces,weights=cces$weight)
load("exit.RData")
w.exit <- svydesign(~1,data=exit,weights=exit$weight)
load("pew.RData")
w.pew <- svydesign(~1,data=pew,weights=pew$weight)

ans.ci <- function(x){
  cis <- confint(x)
  ans <- cbind(as.matrix(x), cis)
  ans <- round(ans, digits=4)
  ans <- as.data.frame(ans)
  names(ans) <- c("est", "lo", "hi")
  ans
}

a.names <- c("Liberal ideology ", "Democratic party ID ", "Presidential vote for Democrat*", 
             "Support marriage equality ", "Support LGB adoption rights*", "Support LGB job protections ",
             "Increase spending and services ", "Reduce defense spending ", "Support government insurance ", "Support guaranteed jobs ", "Support aid to blacks ", "Increase spending on Social Security ", "Increase spending on public schools ",  "Increase spending on science ",  "Increase spending on welfare ", "Increase spending on child care ",   "Increase spending on environment ",  "Oppose death penalty ", "Make it harder to buy guns ", "Increase immigration levels ", "Pro-choice abortion view ", "Support larger government ", "Lower racial resentment ", "Support affirmative action "
)
c.names <- c("Liberal ideology ", "Democratic party ID ", "Presidential vote for Democrat*",
             "Support marriage equality*", "Marriage equality important ", "Oppose transgender military ban*", 
             "Support more gun control ", "Pro-choice abortion views ", "Support environmental protections ", "Oppose repealing Obamacare*", "More liberal immigration views (2016) ", "More liberal immigration views (2018) ", "More liberal racial attitudes ", "Lower racial resentment ", "Support infrastructure spending*", "Support raising minimum wage*", "Oppose Gorsuch confirmation*", "Oppose Kavanaugh confirmation*", "Support Russia sanctions*","Oppose Jerusalem as capital*", "Oppose Keystone pipeline*", "Support Paris climate agreement*", "Oppose TPP withdrawal*", "Support clean power rules*", "Support Iran nuclear deal*", "Oppose travel ban*", "Oppose cutting regulations*"
)
p.names <- c("Liberal ideology ", "Democratic party ID ",
             "Support marriage equality ", "Support LGB adoption rights ",
             "Support larger government*", "Support gun control*", "Believe immigrants strengthen US*", "Benefits don't go far enough*"
)
e.names <- c("Liberal ideology ", "Democratic party ID ", "House vote for Democrat*")





#################################################################################################################
# LGBT DISTINCTIVENESS FROM STRAIGHT CISGENDER AMERICANS
#################################################################################################################

# ANES REGRESSIONS
a.genl.ideo.dum <- svyglm(genl.ideo.c~lgbdummy+age+race+female+married+income+educ+rel.atd+region+year, w.anes)
a.genl.pid7.dum <- svyglm(genl.pid7.c~lgbdummy+age+race+female+married+income+educ+rel.atd+region+year, w.anes)
a.genl.vote.dum <- svyglm(genl.vote.b~lgbdummy+age+race+female+married+income+educ+rel.atd+region+year, w.anes, family=binomial(link="logit"))

a.lgbt.marr.dum <- svyglm(lgbt.marr.c~lgbdummy+age+race+female+married+income+educ+rel.atd+region+year, w.anes)
a.lgbt.adpt.dum <- svyglm(lgbt.adpt.b~lgbdummy+age+race+female+married+income+educ+rel.atd+region+year, w.anes, family=binomial(link="logit"))
a.lgbt.jobs.dum <- svyglm(lgbt.jobs.c~lgbdummy+age+race+female+married+income+educ+rel.atd+region+year, w.anes)

a.atts.serv.dum <- svyglm(atts.serv.c~lgbdummy+age+race+female+married+income+educ+rel.atd+region+year, w.anes)
a.atts.defs.dum <- svyglm(atts.defs.c~lgbdummy+age+race+female+married+income+educ+rel.atd+region+year, w.anes)
a.atts.hlth.dum <- svyglm(atts.hlth.c~lgbdummy+age+race+female+married+income+educ+rel.atd+region+year, w.anes)
a.atts.guar.dum <- svyglm(atts.guar.c~lgbdummy+age+race+female+married+income+educ+rel.atd+region+year, w.anes)
a.atts.aidb.dum <- svyglm(atts.aidb.c~lgbdummy+age+race+female+married+income+educ+rel.atd+region+year, w.anes)
a.atts.socs.dum <- svyglm(atts.socs.c~lgbdummy+age+race+female+married+income+educ+rel.atd+region+year, w.anes)
a.atts.schl.dum <- svyglm(atts.schl.c~lgbdummy+age+race+female+married+income+educ+rel.atd+region+year, w.anes)
a.atts.scie.dum <- svyglm(atts.scie.c~lgbdummy+age+race+female+married+income+educ+rel.atd+region+year, w.anes)
a.atts.wlfr.dum <- svyglm(atts.wlfr.c~lgbdummy+age+race+female+married+income+educ+rel.atd+region+year, w.anes)
a.atts.chil.dum <- svyglm(atts.chil.c~lgbdummy+age+race+female+married+income+educ+rel.atd+region+year, w.anes)
a.atts.envn.dum <- svyglm(atts.envn.c~lgbdummy+age+race+female+married+income+educ+rel.atd+region+year, w.anes)
a.atts.dpen.dum <- svyglm(atts.dpen.c~lgbdummy+age+race+female+married+income+educ+rel.atd+region+year, w.anes)
a.atts.guns.dum <- svyglm(atts.guns.c~lgbdummy+age+race+female+married+income+educ+rel.atd+region+year, w.anes)
a.atts.imlv.dum <- svyglm(atts.imlv.c~lgbdummy+age+race+female+married+income+educ+rel.atd+region+year, w.anes)
a.atts.abrt.dum <- svyglm(atts.abrt.c~lgbdummy+age+race+female+married+income+educ+rel.atd+region+year, w.anes)
a.atts.lgov.dum <- svyglm(atts.lgov.c~lgbdummy+age+race+female+married+income+educ+rel.atd+region+year, w.anes)
a.atts.rsnt.dum <- svyglm(atts.rsnt.c~lgbdummy+age+race+female+married+income+educ+rel.atd+region+year, w.anes)
a.atts.aact.dum <- svyglm(atts.aact.c~lgbdummy+age+race+female+married+income+educ+rel.atd+region+year, w.anes)

# ANES SIMULATIONS
svymean(~age+female+married+income+educ+rel.atd+region+year, w.anes, na.rm=T)
p.anes.0 <- data.frame(lgbdummy=0, age=4.71, race="1White", female=1, married=1, income=3.06, educ=2.91, rel.atd=2.62, region="1South", year="2012") 
p.anes.1 <- data.frame(lgbdummy=1, age=4.71, race="1White", female=1, married=1, income=3.06, educ=2.91, rel.atd= 2.62, region="1South", year="2012") 

a.rgs.dum <- list(a.genl.ideo.dum, a.genl.pid7.dum, a.genl.vote.dum, 
                  a.lgbt.marr.dum, a.lgbt.adpt.dum, a.lgbt.jobs.dum, 
                  a.atts.serv.dum, a.atts.defs.dum, a.atts.hlth.dum, a.atts.guar.dum, a.atts.aidb.dum, a.atts.socs.dum, a.atts.schl.dum, a.atts.scie.dum, a.atts.wlfr.dum, a.atts.chil.dum, a.atts.envn.dum, a.atts.dpen.dum, a.atts.guns.dum, a.atts.imlv.dum, a.atts.abrt.dum, a.atts.lgov.dum, a.atts.rsnt.dum, a.atts.aact.dum)
a.fds.dum <- as.data.frame(matrix(NA, nrow=length(a.rgs.dum), ncol=4))
names(a.fds.dum) <- c("fd", "lo", "hi", "dv")
for(i in 1:length(a.rgs.dum)){
  mod <- a.rgs.dum[[i]]
  fds <- ans.ci(predict(mod, p.anes.1, type="response")-predict(mod, p.anes.0, type="response"))
  a.fds.dum[i,1:3] <- fds
  a.fds.dum[i,4] <- as.character(mod$terms[[2]])
}
a.fds.dum$cat <- c(rep("1genl", 3), rep("2lgbt", 3), rep("3atts", 18))
a.fds.dum$dv <- a.names



# CCES REGRESSIONS
c.genl.ideo.dum <- svyglm(genl.ideo.c~lgbtdummy+age+race+female+married+income+educ+rel.atd+region+year, w.cces)
c.genl.pid7.dum <- svyglm(genl.pid7.c~lgbtdummy+age+race+female+married+income+educ+rel.atd+region+year, w.cces)
c.genl.vote.dum <- svyglm(genl.vote.b~lgbtdummy+age+race+female+married+income+educ+rel.atd+region, w.cces, family=binomial(link="logit"))

c.lgbt.marr.dum <- svyglm(lgbt.marr.b~lgbtdummy+age+race+female+married+income+educ+rel.atd+region, w.cces, family=binomial(link="logit"))
c.lgbt.mimp.dum <- svyglm(lgbt.mimp.c~lgbtdummy+age+race+female+married+income+educ+rel.atd+region, w.cces)
c.lgbt.tban.dum <- svyglm(lgbt.tban.b~lgbtdummy+age+race+female+married+income+educ+rel.atd+region, w.cces, family=binomial(link="logit"))

c.atts.guns.dum <- svyglm(atts.guns.c~lgbtdummy+age+race+female+married+income+educ+rel.atd+region+year, w.cces)
c.atts.abrt.dum <- svyglm(atts.abrt.c~lgbtdummy+age+race+female+married+income+educ+rel.atd+region+year, w.cces)
c.atts.envn.dum <- svyglm(atts.envn.c~lgbtdummy+age+race+female+married+income+educ+rel.atd+region+year, w.cces)
c.atts.ocre.dum <- svyglm(atts.ocre.b~lgbtdummy+age+race+female+married+income+educ+rel.atd+region, w.cces, family=binomial(link="logit"))
c.atts.im16.dum <- svyglm(atts.im16.c~lgbtdummy+age+race+female+married+income+educ+rel.atd+region, w.cces)
c.atts.im18.dum <- svyglm(atts.im18.c~lgbtdummy+age+race+female+married+income+educ+rel.atd+region, w.cces)
c.atts.racm.dum <- svyglm(atts.racm.c~lgbtdummy+age+race+female+married+income+educ+rel.atd+region, w.cces)
c.atts.rsnt.dum <- svyglm(atts.rsnt.c~lgbtdummy+age+race+female+married+income+educ+rel.atd+region, w.cces)
c.atts.infr.dum <- svyglm(atts.infr.b~lgbtdummy+age+race+female+married+income+educ+rel.atd+region, w.cces, family=binomial(link="logit"))
c.atts.mwag.dum <- svyglm(atts.mwag.b~lgbtdummy+age+race+female+married+income+educ+rel.atd+region, w.cces, family=binomial(link="logit"))
c.atts.gors.dum <- svyglm(atts.gors.b~lgbtdummy+age+race+female+married+income+educ+rel.atd+region, w.cces, family=binomial(link="logit"))
c.atts.kavn.dum <- svyglm(atts.kavn.b~lgbtdummy+age+race+female+married+income+educ+rel.atd+region, w.cces, family=binomial(link="logit"))
c.atts.russ.dum <- svyglm(atts.russ.b~lgbtdummy+age+race+female+married+income+educ+rel.atd+region, w.cces, family=binomial(link="logit"))
c.atts.isrl.dum <- svyglm(atts.isrl.b~lgbtdummy+age+race+female+married+income+educ+rel.atd+region, w.cces, family=binomial(link="logit"))
c.atts.keys.dum <- svyglm(atts.keys.b~lgbtdummy+age+race+female+married+income+educ+rel.atd+region, w.cces, family=binomial(link="logit"))
c.atts.clim.dum <- svyglm(atts.clim.b~lgbtdummy+age+race+female+married+income+educ+rel.atd+region, w.cces, family=binomial(link="logit"))
c.atts.tp18.dum <- svyglm(atts.tp18.b~lgbtdummy+age+race+female+married+income+educ+rel.atd+region, w.cces, family=binomial(link="logit"))
c.atts.powr.dum <- svyglm(atts.powr.b~lgbtdummy+age+race+female+married+income+educ+rel.atd+region, w.cces, family=binomial(link="logit"))
c.atts.ir18.dum <- svyglm(atts.ir18.b~lgbtdummy+age+race+female+married+income+educ+rel.atd+region, w.cces, family=binomial(link="logit"))
c.atts.mban.dum <- svyglm(atts.mban.b~lgbtdummy+age+race+female+married+income+educ+rel.atd+region, w.cces, family=binomial(link="logit"))
c.atts.regs.dum <- svyglm(atts.regs.b~lgbtdummy+age+race+female+married+income+educ+rel.atd+region, w.cces, family=binomial(link="logit"))

# CCES SIMULATIONS
svymean(~age+female+married+income+educ+rel.atd+region, w.cces, na.rm=T)
p.cces.0 <- data.frame(lgbtdummy=0, age=4.72, race="1White", female=1, married=1, income=2.75, educ=2.90, rel.atd=2.91, region="1South", year="2016") 
p.cces.1 <- data.frame(lgbtdummy=1, age=4.72, race="1White", female=1, married=1, income=2.75, educ=2.90, rel.atd= 2.91, region="1South", year="2016") 

c.rgs.dum <- list(c.genl.ideo.dum, c.genl.pid7.dum, c.genl.vote.dum, 
                  c.lgbt.marr.dum, c.lgbt.mimp.dum, c.lgbt.tban.dum, 
                  c.atts.guns.dum, c.atts.abrt.dum, c.atts.envn.dum, c.atts.ocre.dum, c.atts.im16.dum, c.atts.im18.dum, c.atts.racm.dum, 
                  c.atts.rsnt.dum, c.atts.infr.dum, c.atts.mwag.dum, c.atts.gors.dum, c.atts.kavn.dum, 
                  c.atts.russ.dum, c.atts.isrl.dum, c.atts.keys.dum, c.atts.clim.dum, c.atts.tp18.dum, c.atts.powr.dum, 
                  c.atts.ir18.dum, c.atts.mban.dum, c.atts.regs.dum)
c.fds.dum <- as.data.frame(matrix(NA, nrow=length(c.rgs.dum), ncol=4))
names(c.fds.dum) <- c("fd", "lo", "hi", "dv")
for(i in 1:length(c.rgs.dum)){
  mod <- c.rgs.dum[[i]]
  fds <- ans.ci(predict(mod, p.cces.1, type="response")-predict(mod, p.cces.0, type="response"))
  c.fds.dum[i,1:3] <- fds
  c.fds.dum[i,4] <- as.character(mod$terms[[2]])
}
c.fds.dum$cat <- c(rep("1genl", 3), rep("2lgbt", 3), rep("3atts", 21))
c.fds.dum$dv <- c.names



# AP VOTECAST EXIT POLL REGRESSIONS
e.genl.ideo.dum <- svyglm(genl.ideo.c~lgbtdummy+age+race+female+married+income+educ+rel.atd+region, w.exit)
e.genl.pid5.dum <- svyglm(genl.pid5.c~lgbtdummy+age+race+female+married+income+educ+rel.atd+region, w.exit)
e.genl.voth.dum <- svyglm(genl.voth.b~lgbtdummy+age+race+female+married+income+educ+rel.atd+region, w.exit, family=binomial(link="logit"))
e.rgs.dum <- list(e.genl.ideo.dum, e.genl.pid5.dum, e.genl.voth.dum)

# AP VOTECAST EXIT POLL SIMULATIONS
svymean(~age+female+married+income+educ+rel.atd+region, w.exit, na.rm=T)
p.exit.0 <- data.frame(lgbtdummy=0, age=3.88, race="1White", female=1, married=1, income=3.00, educ=2.09, rel.atd=2.55, region="1South") 
p.exit.1 <- data.frame(lgbtdummy=1, age=3.88, race="1White", female=1, married=1, income=3.00, educ=2.09, rel.atd=2.55, region="1South") 

e.fds.dum <- as.data.frame(matrix(NA, nrow=length(e.rgs.dum), ncol=4))
names(e.fds.dum) <- c("fd", "lo", "hi", "dv")
for(i in 1:length(e.rgs.dum)){
  mod <- e.rgs.dum[[i]]
  fds <- ans.ci(predict(mod, p.exit.1, type="response")-predict(mod, p.exit.0, type="response"))
  e.fds.dum[i,1:3] <- fds
  e.fds.dum[i,4] <- as.character(mod$terms[[2]])
}
e.fds.dum$cat <- "1genl"
e.fds.dum$dv <- e.names





#################################################################################################################
# DIVERSITY WITHIN THE LGBT COMMUNITY
#################################################################################################################

# ANES REGRESSIONS
a.genl.ideo.div <- svyglm(genl.ideo.c~lgb2+age+race+female+married+income+educ+rel.atd+region+year, w.anes)
a.genl.pid7.div <- svyglm(genl.pid7.c~lgb2+age+race+female+married+income+educ+rel.atd+region+year, w.anes)
a.genl.vote.div <- svyglm(genl.vote.b~lgb2+age+race+female+married+income+educ+rel.atd+region+year, w.anes, family=binomial(link="logit"))

a.lgbt.marr.div <- svyglm(lgbt.marr.c~lgb2+age+race+female+married+income+educ+rel.atd+region+year, w.anes)
a.lgbt.adpt.div <- svyglm(lgbt.adpt.b~lgb2+age+race+female+married+income+educ+rel.atd+region+year, w.anes, family=binomial(link="logit"))
a.lgbt.jobs.div <- svyglm(lgbt.jobs.c~lgb2+age+race+female+married+income+educ+rel.atd+region+year, w.anes)

a.atts.serv.div <- svyglm(atts.serv.c~lgb2+age+race+female+married+income+educ+rel.atd+region+year, w.anes)
a.atts.defs.div <- svyglm(atts.defs.c~lgb2+age+race+female+married+income+educ+rel.atd+region+year, w.anes)
a.atts.hlth.div <- svyglm(atts.hlth.c~lgb2+age+race+female+married+income+educ+rel.atd+region+year, w.anes)
a.atts.guar.div <- svyglm(atts.guar.c~lgb2+age+race+female+married+income+educ+rel.atd+region+year, w.anes)
a.atts.aidb.div <- svyglm(atts.aidb.c~lgb2+age+race+female+married+income+educ+rel.atd+region+year, w.anes)
a.atts.socs.div <- svyglm(atts.socs.c~lgb2+age+race+female+married+income+educ+rel.atd+region+year, w.anes)
a.atts.schl.div <- svyglm(atts.schl.c~lgb2+age+race+female+married+income+educ+rel.atd+region+year, w.anes)
a.atts.scie.div <- svyglm(atts.scie.c~lgb2+age+race+female+married+income+educ+rel.atd+region+year, w.anes)
a.atts.wlfr.div <- svyglm(atts.wlfr.c~lgb2+age+race+female+married+income+educ+rel.atd+region+year, w.anes)
a.atts.chil.div <- svyglm(atts.chil.c~lgb2+age+race+female+married+income+educ+rel.atd+region+year, w.anes)
a.atts.envn.div <- svyglm(atts.envn.c~lgb2+age+race+female+married+income+educ+rel.atd+region+year, w.anes)
a.atts.dpen.div <- svyglm(atts.dpen.c~lgb2+age+race+female+married+income+educ+rel.atd+region+year, w.anes)
a.atts.guns.div <- svyglm(atts.guns.c~lgb2+age+race+female+married+income+educ+rel.atd+region+year, w.anes)
a.atts.imlv.div <- svyglm(atts.imlv.c~lgb2+age+race+female+married+income+educ+rel.atd+region+year, w.anes)
a.atts.abrt.div <- svyglm(atts.abrt.c~lgb2+age+race+female+married+income+educ+rel.atd+region+year, w.anes)
a.atts.lgov.div <- svyglm(atts.lgov.c~lgb2+age+race+female+married+income+educ+rel.atd+region+year, w.anes)
a.atts.rsnt.div <- svyglm(atts.rsnt.c~lgb2+age+race+female+married+income+educ+rel.atd+region+year, w.anes)
a.atts.aact.div <- svyglm(atts.aact.c~lgb2+age+race+female+married+income+educ+rel.atd+region+year, w.anes)


# ANES SIMULATIONS
p.anes.lg <- data.frame(lgb2="2Gay", age=4.71, race="1White", female=1, married=1, income=3.06, educ=2.91, rel.atd=2.62, region="1South", year="2012") 
p.anes.bi <- data.frame(lgb2="3Bi", age=4.71, race="1White", female=1, married=1, income=3.06, educ=2.91, rel.atd= 2.62, region="1South", year="2012") 

a.rgs.div <- list(a.genl.ideo.div, a.genl.pid7.div, a.genl.vote.div, 
                  a.lgbt.marr.div, a.lgbt.adpt.div, a.lgbt.jobs.div,
                  a.atts.serv.div, a.atts.defs.div, a.atts.hlth.div, a.atts.guar.div, a.atts.aidb.div, a.atts.socs.div, a.atts.schl.div, 		
                  a.atts.scie.div, a.atts.wlfr.div, a.atts.chil.div, a.atts.envn.div, a.atts.dpen.div, a.atts.guns.div, a.atts.imlv.div, 
                  a.atts.abrt.div, a.atts.lgov.div, a.atts.rsnt.div, a.atts.aact.div)
a.fds.div <- as.data.frame(matrix(NA, nrow=length(a.rgs.div), ncol=4))
names(a.fds.div) <- c("fd", "lo", "hi", "dv")

for(i in 1:length(a.rgs.div)){
  mod <- a.rgs.div[[i]]
  fds <- ans.ci(predict(mod, p.anes.bi, type="response")-predict(mod, p.anes.lg, type="response"))
  a.fds.div[i,1:3] <- fds
  a.fds.div[i,4] <- as.character(mod$terms[[2]])
}
a.fds.div$cat <- c(rep("1genl", 3), rep("2lgbt", 3), rep("3atts", 18))
a.fds.div$dv <- a.names



# CCES REGRESSIONS
c.genl.ideo.div <- svyglm(genl.ideo.c~lgb2+trans+age+race+female+married+income+educ+rel.atd+region+year, w.cces)
c.genl.pid7.div <- svyglm(genl.pid7.c~lgb2+trans+age+race+female+married+income+educ+rel.atd+region+year, w.cces)
c.genl.vote.div <- svyglm(genl.vote.b~lgb2+trans+age+race+female+married+income+educ+rel.atd+region, w.cces, family=binomial(link="logit"))

c.lgbt.marr.div <- svyglm(lgbt.marr.b~lgb2+trans+age+race+female+married+income+educ+rel.atd+region, w.cces, family=binomial(link="logit"))
c.lgbt.mimp.div <- svyglm(lgbt.mimp.c~lgb2+trans+age+race+female+married+income+educ+rel.atd+region, w.cces)
c.lgbt.tban.div <- svyglm(lgbt.tban.b~lgb2+trans+age+race+female+married+income+educ+rel.atd+region, w.cces, family=binomial(link="logit"))

c.atts.guns.div <- svyglm(atts.guns.c~lgb2+trans+age+race+female+married+income+educ+rel.atd+region+year, w.cces)
c.atts.abrt.div <- svyglm(atts.abrt.c~lgb2+trans+age+race+female+married+income+educ+rel.atd+region+year, w.cces)
c.atts.envn.div <- svyglm(atts.envn.c~lgb2+trans+age+race+female+married+income+educ+rel.atd+region+year, w.cces)
c.atts.ocre.div <- svyglm(atts.ocre.b~lgb2+trans+age+race+female+married+income+educ+rel.atd+region, w.cces, family=binomial(link="logit"))
c.atts.im16.div <- svyglm(atts.im16.c~lgb2+trans+age+race+female+married+income+educ+rel.atd+region, w.cces)
c.atts.im18.div <- svyglm(atts.im18.c~lgb2+trans+age+race+female+married+income+educ+rel.atd+region, w.cces)
c.atts.racm.div <- svyglm(atts.racm.c~lgb2+trans+age+race+female+married+income+educ+rel.atd+region, w.cces)
c.atts.rsnt.div <- svyglm(atts.rsnt.c~lgb2+trans+age+race+female+married+income+educ+rel.atd+region, w.cces)
c.atts.infr.div <- svyglm(atts.infr.b~lgb2+trans+age+race+female+married+income+educ+rel.atd+region, w.cces, family=binomial(link="logit"))
c.atts.mwag.div <- svyglm(atts.mwag.b~lgb2+trans+age+race+female+married+income+educ+rel.atd+region, w.cces, family=binomial(link="logit"))
c.atts.gors.div <- svyglm(atts.gors.b~lgb2+trans+age+race+female+married+income+educ+rel.atd+region, w.cces, family=binomial(link="logit"))
c.atts.kavn.div <- svyglm(atts.kavn.b~lgb2+trans+age+race+female+married+income+educ+rel.atd+region, w.cces, family=binomial(link="logit"))
c.atts.russ.div <- svyglm(atts.russ.b~lgb2+trans+age+race+female+married+income+educ+rel.atd+region, w.cces, family=binomial(link="logit"))
c.atts.isrl.div <- svyglm(atts.isrl.b~lgb2+trans+age+race+female+married+income+educ+rel.atd+region, w.cces, family=binomial(link="logit"))
c.atts.keys.div <- svyglm(atts.keys.b~lgb2+trans+age+race+female+married+income+educ+rel.atd+region, w.cces, family=binomial(link="logit"))
c.atts.clim.div <- svyglm(atts.clim.b~lgb2+trans+age+race+female+married+income+educ+rel.atd+region, w.cces, family=binomial(link="logit"))
c.atts.tp18.div <- svyglm(atts.tp18.b~lgb2+trans+age+race+female+married+income+educ+rel.atd+region, w.cces, family=binomial(link="logit"))
c.atts.powr.div <- svyglm(atts.powr.b~lgb2+trans+age+race+female+married+income+educ+rel.atd+region, w.cces, family=binomial(link="logit"))
c.atts.ir18.div <- svyglm(atts.ir18.b~lgb2+trans+age+race+female+married+income+educ+rel.atd+region, w.cces, family=binomial(link="logit"))
c.atts.mban.div <- svyglm(atts.mban.b~lgb2+trans+age+race+female+married+income+educ+rel.atd+region, w.cces, family=binomial(link="logit"))
c.atts.regs.div <- svyglm(atts.regs.b~lgb2+trans+age+race+female+married+income+educ+rel.atd+region, w.cces, family=binomial(link="logit"))

# CCES SIMULATIONS
svymean(~age+female+married+income+educ+rel.atd+region, w.cces, na.rm=T)
p.cces.lg <- data.frame(lgb2="2Gay", trans=0, age=4.72, race="1White", female=1, married=1, income=2.75, educ=2.90, rel.atd=2.91, region="1South", year="2016") 
p.cces.bi <- data.frame(lgb2="3Bi", trans=0, age=4.72, race="1White", female=1, married=1, income=2.75, educ=2.90, rel.atd= 2.91, region="1South", year="2016") 
p.cces.tr <- data.frame(lgb2="1Straight", trans=1, age=4.72, race="1White", female=1, married=1, income=2.75, educ=2.90, rel.atd=2.91, region="1South", year="2016") 


c.rgs.div <- list(c.genl.ideo.div, c.genl.pid7.div, c.genl.vote.div, 
                  c.lgbt.marr.div, c.lgbt.mimp.div, c.lgbt.tban.div, 
                  c.atts.guns.div, c.atts.abrt.div, c.atts.envn.div, c.atts.ocre.div, c.atts.im16.div, c.atts.im18.div, c.atts.racm.div, 
                  c.atts.rsnt.div, c.atts.infr.div, c.atts.mwag.div, c.atts.gors.div, c.atts.kavn.div, 
                  c.atts.russ.div, c.atts.isrl.div, c.atts.keys.div, c.atts.clim.div, c.atts.tp18.div, c.atts.powr.div, 
                  c.atts.ir18.div, c.atts.mban.div, c.atts.regs.div)
c.fds.div <- as.data.frame(matrix(NA, nrow=length(c.rgs.div), ncol=4))
names(c.fds.div) <- c("fd", "lo", "hi", "dv")

c.fds.div.bi <- c.fds.div
c.fds.div.tr <- c.fds.div


for(i in 1:length(c.rgs.div)){
  mod <- c.rgs.div[[i]]
  fds <- ans.ci(predict(mod, p.cces.bi, type="response")-predict(mod, p.cces.lg, type="response"))
  c.fds.div.bi[i,1:3] <- fds
  c.fds.div.bi[i,4] <- as.character(mod$terms[[2]])
  
  fds <- ans.ci(predict(mod, p.cces.tr, type="response")-predict(mod, p.cces.lg, type="response"))
  c.fds.div.tr[i,1:3] <- fds
  c.fds.div.tr[i,4] <- as.character(mod$terms[[2]])
}
c.fds.div.bi$cat <- c(rep("1genl", 3), rep("2lgbt", 3), rep("3atts", 21))
c.fds.div.tr$cat <- c(rep("1genl", 3), rep("2lgbt", 3), rep("3atts", 21))
c.fds.div.bi$dv <- c.names
c.fds.div.tr$dv <- c.names



# PEW REGRESSIONS
p.genl.ideo.div <- svyglm(genl.ideo.c~lgbt+age+race+female+married+income+educ+rel.atd, w.pew)
p.genl.pid3.div <- svyglm(genl.pid3.c~lgbt+age+race+female+married+income+educ+rel.atd, w.pew)

p.lgbt.marr.div <- svyglm(lgbt.marr.c~lgbt+age+race+female+married+income+educ+rel.atd, w.pew)
p.lgbt.adpt.div <- svyglm(lgbt.adpt.c~lgbt+age+race+female+married+income+educ+rel.atd, w.pew)

p.atts.govt.div <- svyglm(atts.govt.b~lgbt+age+race+female+married+income+educ+rel.atd, w.pew, family=binomial(link="logit"))
p.atts.guns.div <- svyglm(atts.guns.b~lgbt+age+race+female+married+income+educ+rel.atd, w.pew, family=binomial(link="logit"))
p.atts.immg.div <- svyglm(atts.immg.b~lgbt+age+race+female+married+income+educ+rel.atd, w.pew, family=binomial(link="logit"))
p.atts.bens.div <- svyglm(atts.bens.b~lgbt+age+race+female+married+income+educ+rel.atd, w.pew, family=binomial(link="logit"))

# PEW SIMULATIONS
svymean(~lgbt+age+female+married+income+educ+rel.atd, w.pew, na.rm=T)
p.pew.lg <- data.frame(lgbt="2Gay", age=3.19, race="1White", female=1, married=0, income=2.40, educ=2.93, rel.atd=2.30)
p.pew.bi <- data.frame(lgbt="3Bi", age= 3.19, race="1White", female=1, married=0, income=2.40, educ=2.93, rel.atd= 2.30)
p.pew.tr <- data.frame(lgbt="4Trans", age=3.19, race="1White", female=1, married=0, income=2.40, educ=2.93, rel.atd=2.30)

p.rgs.div <- list(p.genl.ideo.div, p.genl.pid3.div, 
                  p.lgbt.marr.div, p.lgbt.adpt.div, 
                  p.atts.govt.div, p.atts.guns.div, p.atts.immg.div, p.atts.bens.div)
p.fds.div <- as.data.frame(matrix(NA, nrow=length(p.rgs.div), ncol=4))
names(p.fds.div) <- c("fd", "lo", "hi", "dv")

p.fds.div.bi <- p.fds.div
p.fds.div.tr <- p.fds.div

for(i in 1:length(p.rgs.div)){
  mod <- p.rgs.div[[i]]
  fds <- ans.ci(predict(mod, p.pew.bi, type="response")-predict(mod, p.pew.lg, type="response"))
  p.fds.div.bi[i,1:3] <- fds
  p.fds.div.bi[i,4] <- as.character(mod$terms[[2]])
  
  fds <- ans.ci(predict(mod, p.pew.tr, type="response")-predict(mod, p.pew.lg, type="response"))
  p.fds.div.tr[i,1:3] <- fds
  p.fds.div.tr[i,4] <- as.character(mod$terms[[2]])
}

p.fds.div.bi$cat <- c(rep("1genl", 2), rep("2lgbt", 2), rep("3atts", 4))
p.fds.div.tr$cat <- c(rep("1genl", 2), rep("2lgbt", 2), rep("3atts", 4))
p.fds.div.bi$dv <- p.names
p.fds.div.tr$dv <- p.names




# EXIT REGRESSIONS
e.genl.ideo.div <- svyglm(genl.ideo.c~lgb+trans+age+race+female+married+income+educ+rel.atd+region, w.exit)
e.genl.pid5.div <- svyglm(genl.pid5.c~lgb+trans+age+race+female+married+income+educ+rel.atd+region, w.exit)
e.genl.voth.div <- svyglm(genl.voth.b~lgb+trans+age+race+female+married+income+educ+rel.atd+region, w.exit, family=binomial(link="logit"))

e.rgs.div <- list(e.genl.ideo.div, e.genl.pid5.div, e.genl.voth.div)

svymean(~age+female+married+income+educ+rel.atd+region, w.exit, na.rm=T)
p.exit.0 <- data.frame(lgb=1, trans=0, age=3.88, race="1White", female=1, married=1, income=3.00, educ=2.09, rel.atd=2.55, region="1South") 
p.exit.1 <- data.frame(lgb=0, trans=1, age=3.88, race="1White", female=1, married=1, income=3.00, educ=2.09, rel.atd=2.55, region="1South") 


# EXIT SIMULATIONS
e.fds.div <- as.data.frame(matrix(NA, nrow=length(e.rgs.div), ncol=4))
names(e.fds.div) <- c("fd", "lo", "hi", "dv")
for(i in 1:length(e.rgs.div)){
  mod <- e.rgs.div[[i]]
  fds <- ans.ci(predict(mod, p.exit.1, type="response")-predict(mod, p.exit.0, type="response"))
  e.fds.div[i,1:3] <- fds
  e.fds.div[i,4] <- as.character(mod$terms[[2]])
}

e.fds.div$cat <- "1genl"
e.fds.div$dv <- e.names




#################################################################################################################
# DIFFERENCES BASED ON GENDER [CCES ONLY]
#################################################################################################################

c.genl.ideo.gen <- svyglm(genl.ideo.c~lgb2*female+trans*female+age+race+married+income+educ+rel.atd+region+year, w.cces)
c.genl.pid7.gen <- svyglm(genl.pid7.c~lgb2*female+trans*female+age+race+married+income+educ+rel.atd+region+year, w.cces)
c.genl.vote.gen <- svyglm(genl.vote.b~lgb2*female+trans*female+age+race+married+income+educ+rel.atd+region, w.cces, family=binomial(link="logit"))

c.lgbt.marr.gen <- svyglm(lgbt.marr.b~lgb2*female+trans*female+age+race+married+income+educ+rel.atd+region, w.cces, family=binomial(link="logit"))
c.lgbt.mimp.gen <- svyglm(lgbt.mimp.c~lgb2*female+trans*female+age+race+married+income+educ+rel.atd+region, w.cces)
c.lgbt.tban.gen <- svyglm(lgbt.tban.b~lgb2*female+trans*female+age+race+married+income+educ+rel.atd+region, w.cces, family=binomial(link="logit"))

c.atts.guns.gen <- svyglm(atts.guns.c~lgb2*female+trans*female+age+race+married+income+educ+rel.atd+region+year, w.cces)
c.atts.abrt.gen <- svyglm(atts.abrt.c~lgb2*female+trans*female+age+race+married+income+educ+rel.atd+region+year, w.cces)
c.atts.envn.gen <- svyglm(atts.envn.c~lgb2*female+trans*female+age+race+married+income+educ+rel.atd+region+year, w.cces)
c.atts.ocre.gen <- svyglm(atts.ocre.b~lgb2*female+trans*female+age+race+married+income+educ+rel.atd+region, w.cces, family=binomial(link="logit"))
c.atts.im16.gen <- svyglm(atts.im16.c~lgb2*female+trans*female+age+race+married+income+educ+rel.atd+region, w.cces)
c.atts.im18.gen <- svyglm(atts.im18.c~lgb2*female+trans*female+age+race+married+income+educ+rel.atd+region, w.cces)
c.atts.racm.gen <- svyglm(atts.racm.c~lgb2*female+trans*female+age+race+married+income+educ+rel.atd+region, w.cces)
c.atts.rsnt.gen <- svyglm(atts.rsnt.c~lgb2*female+trans*female+age+race+married+income+educ+rel.atd+region, w.cces)
c.atts.infr.gen <- svyglm(atts.infr.b~lgb2*female+trans*female+age+race+married+income+educ+rel.atd+region, w.cces, family=binomial(link="logit"))
c.atts.mwag.gen <- svyglm(atts.mwag.b~lgb2*female+trans*female+age+race+married+income+educ+rel.atd+region, w.cces, family=binomial(link="logit"))
c.atts.gors.gen <- svyglm(atts.gors.b~lgb2*female+trans*female+age+race+married+income+educ+rel.atd+region, w.cces, family=binomial(link="logit"))
c.atts.kavn.gen <- svyglm(atts.kavn.b~lgb2*female+trans*female+age+race+married+income+educ+rel.atd+region, w.cces, family=binomial(link="logit"))
c.atts.russ.gen <- svyglm(atts.russ.b~lgb2*female+trans*female+age+race+married+income+educ+rel.atd+region, w.cces, family=binomial(link="logit"))
c.atts.isrl.gen <- svyglm(atts.isrl.b~lgb2*female+trans*female+age+race+married+income+educ+rel.atd+region, w.cces, family=binomial(link="logit"))
c.atts.keys.gen <- svyglm(atts.keys.b~lgb2*female+trans*female+age+race+married+income+educ+rel.atd+region, w.cces, family=binomial(link="logit"))
c.atts.clim.gen <- svyglm(atts.clim.b~lgb2*female+trans*female+age+race+married+income+educ+rel.atd+region, w.cces, family=binomial(link="logit"))
c.atts.tp18.gen <- svyglm(atts.tp18.b~lgb2*female+trans*female+age+race+married+income+educ+rel.atd+region, w.cces, family=binomial(link="logit"))
c.atts.powr.gen <- svyglm(atts.powr.b~lgb2*female+trans*female+age+race+married+income+educ+rel.atd+region, w.cces, family=binomial(link="logit"))
c.atts.ir18.gen <- svyglm(atts.ir18.b~lgb2*female+trans*female+age+race+married+income+educ+rel.atd+region, w.cces, family=binomial(link="logit"))
c.atts.mban.gen <- svyglm(atts.mban.b~lgb2*female+trans*female+age+race+married+income+educ+rel.atd+region, w.cces, family=binomial(link="logit"))
c.atts.regs.gen <- svyglm(atts.regs.b~lgb2*female+trans*female+age+race+married+income+educ+rel.atd+region, w.cces, family=binomial(link="logit"))


# simulations
svymean(~age+female+married+income+educ+rel.atd+region, w.cces, na.rm=T)
p.cces.lg.m <- data.frame(lgb2="2Gay", trans=0, age=4.72, race="1White", female=0, married=1, income=2.75, educ=2.90, rel.atd=2.91, region="1South", year="2016") 
p.cces.lg.f <- data.frame(lgb2="2Gay", trans=0, age=4.72, race="1White", female=1, married=1, income=2.75, educ=2.90, rel.atd=2.91, region="1South", year="2016") 
p.cces.bi.m <- data.frame(lgb2="3Bi", trans=0, age=4.72, race="1White", female=0, married=1, income=2.75, educ=2.90, rel.atd= 2.91, region="1South", year="2016") 
p.cces.bi.f <- data.frame(lgb2="3Bi", trans=0, age=4.72, race="1White", female=1, married=1, income=2.75, educ=2.90, rel.atd= 2.91, region="1South", year="2016") 
p.cces.tr.m <- data.frame(lgb2="1Straight", trans=1, age=4.72, race="1White", female=0, married=1, income=2.75, educ=2.90, rel.atd=2.91, region="1South", year="2016") 
p.cces.tr.f <- data.frame(lgb2="1Straight", trans=1, age=4.72, race="1White", female=1, married=1, income=2.75, educ=2.90, rel.atd=2.91, region="1South", year="2016") 


c.rgs.gen <- list(c.genl.ideo.gen, c.genl.pid7.gen, c.genl.vote.gen, 
                  c.lgbt.marr.gen, c.lgbt.mimp.gen, c.lgbt.tban.gen, 
                  c.atts.guns.gen, c.atts.abrt.gen, c.atts.envn.gen, c.atts.ocre.gen, c.atts.im16.gen, 								c.atts.im18.gen, c.atts.racm.gen, 
                  c.atts.rsnt.gen, c.atts.infr.gen, c.atts.mwag.gen, c.atts.gors.gen, c.atts.kavn.gen, 
                  c.atts.russ.gen, c.atts.isrl.gen, c.atts.keys.gen, c.atts.clim.gen, c.atts.tp18.gen, 			
                  c.atts.powr.gen, 
                  c.atts.ir18.gen, c.atts.mban.gen, c.atts.regs.gen)
c.evs.gen <- as.data.frame(matrix(NA, nrow=length(c.rgs.gen), ncol=6))
names(c.evs.gen) <- c("ev", "lo", "hi", "dv", "gender", "lgbt")
c.fds.gen <- as.data.frame(matrix(NA, nrow=length(c.rgs.gen), ncol=5))
names(c.fds.gen) <- c("fd", "lo", "hi", "dv", "lgbt")


c.evs.gen.lg.m <- c.evs.gen
c.evs.gen.lg.f <- c.evs.gen
c.evs.gen.bi.m <- c.evs.gen
c.evs.gen.bi.f <- c.evs.gen
c.evs.gen.tr.m <- c.evs.gen
c.evs.gen.tr.f <- c.evs.gen

c.fds.gen.lg <- c.fds.gen
c.fds.gen.bi <- c.fds.gen
c.fds.gen.tr <- c.fds.gen



for(i in 1:length(c.rgs.gen)){
  mod <- c.rgs.gen[[i]]
  
  c.evs.gen.lg.m[i,1:3] <- ans.ci(predict(mod, p.cces.lg.m, type="response"))
  c.evs.gen.lg.m[i,4] <- as.character(mod$terms[[2]])
  c.evs.gen.lg.m[i,5:6] <- c("Men", "LG")
  c.evs.gen.lg.f[i,1:3] <- ans.ci(predict(mod, p.cces.lg.f, type="response"))
  c.evs.gen.lg.f[i,4] <- as.character(mod$terms[[2]])
  c.evs.gen.lg.f[i,5:6] <- c("Women", "LG")
  
  c.fds.gen.lg[i,1:3] <- ans.ci(predict(mod, p.cces.lg.f, type="response")-predict(mod, p.cces.lg.m, type="response"))
  c.fds.gen.lg[i,4] <- as.character(mod$terms[[2]])
  c.fds.gen.lg[i,5] <- "LG"
  
  
  c.evs.gen.bi.m[i,1:3] <- ans.ci(predict(mod, p.cces.bi.m, type="response"))
  c.evs.gen.bi.m[i,4] <- as.character(mod$terms[[2]])
  c.evs.gen.bi.m[i,5:6] <- c("Men", "Bi")
  c.evs.gen.bi.f[i,1:3] <- ans.ci(predict(mod, p.cces.bi.f, type="response"))
  c.evs.gen.bi.f[i,4] <- as.character(mod$terms[[2]])
  c.evs.gen.bi.f[i,5:6] <- c("Women", "Bi")
  
  c.fds.gen.bi[i,1:3] <- ans.ci(predict(mod, p.cces.bi.f, type="response")-predict(mod, p.cces.bi.m, type="response"))
  c.fds.gen.bi[i,4] <- as.character(mod$terms[[2]])
  c.fds.gen.bi[i,5] <- "Bi"
  
  
  c.evs.gen.tr.m[i,1:3] <- ans.ci(predict(mod, p.cces.tr.m, type="response"))
  c.evs.gen.tr.m[i,4] <- as.character(mod$terms[[2]])
  c.evs.gen.tr.m[i,5:6] <- c("Men", "Tr")
  c.evs.gen.tr.f[i,1:3] <- ans.ci(predict(mod, p.cces.tr.f, type="response"))
  c.evs.gen.tr.f[i,4] <- as.character(mod$terms[[2]])
  c.evs.gen.tr.f[i,5:6] <- c("Women", "Tr")
  
  c.fds.gen.tr[i,1:3] <- ans.ci(predict(mod, p.cces.tr.f, type="response")-predict(mod, p.cces.tr.m, type="response"))
  c.fds.gen.tr[i,4] <- as.character(mod$terms[[2]])
  c.fds.gen.tr[i,5] <- "Tr"
}

averages <- matrix(NA, ncol=1, nrow=27)
for(i in 1:27){
  ave <- mean(c.evs.gen.lg.m$ev[i], c.evs.gen.lg.f$ev[i],  c.evs.gen.bi.m$ev[i], c.evs.gen.bi.f$ev[i],  c.evs.gen.tr.m$ev[i], c.evs.gen.tr.f$ev[i])
  averages[i] <- ave
}
order <- c(1:3, 3,1,2, rank(averages[7:27], ties.method="random"))

c.evs.gen <- rbind(c.evs.gen.lg.f, c.evs.gen.lg.m, c.evs.gen.bi.f, c.evs.gen.bi.m, c.evs.gen.tr.f, c.evs.gen.tr.m)
c.evs.gen$dv <- c.names
c.evs.gen$cat <- c(rep("1genl", 3), rep("2lgbt", 3), rep("3atts", 21))
c.evs.gen$order <- order
rm(averages, order)
c.evs.gen$gender <- as.factor(c.evs.gen$gender)